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Abstract. We study the inverse electrostatic and elasticity prob- 
lems associated with Poisson and Navier equations. The unique- 
ness of solutions of these problems is proved for piecewise constant 
electric charge and internal stress distributions having a checkered 
structure: they are constant on rectangular blocks. Such distribu- 
tions appear naturally in practical applications. We also discuss 
computational challenges arising in the numerical implementation 
of our method. 



I. Introduction and main results 

1.1. Direct electrostatic and elasticity problems. The Poisson 
and Navier (also known as Lame) partial differential equations of el- 
liptic type are commonly used for direct problems in electrostatics and 
elasticity theory Let Q be a bounded domain in IR n , n = 2,3,...; 
the cases n = 2, 3 are the most physically interesting. It is assumed 
throughout the paper that Q has a piecewise smooth boundary F" = dQ. 

In the direct formulation of the electrostatic problem, the Poisson 
equation for a real valued function u(x), called the electric potential 
distribution, 

Au = f(x), (1.1.1) 
is solved in the domain Q with a known distribution of the electric 
charge density, —f(x), and with definite boundary conditions set on 
r. The boundary conditions may be formulated either in the form of 
potential values (Dirichlet conditions) 

w|r = 0i (1-1.2) 
or in terms of the electric field (Neumann conditions) , 

(V«,i/)| r = 02. (1-1.3) 

Here v = (vi, . . . , v n ) is the unit outer normal vector to V and (-, •) is the 
Euclidean scalar product in W 1 . Different parts of Y may have different 
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types of boundary conditions, and at any part of T only one boundary 
condition may be set (which may be a linear combination of Dirichlet 
and Neumann conditions), so that the problem is not overconstrained. 

The direct formulation of the elasticity problem is described by the 
Navier equation 

AU + agmddivU = F{x), (1.1.4) 

for a vector- valued function U : Q — > W 1 , called the displacement field. 
Here F(x) = — 2(1 +^ J-~(x), where T is the distribution of body forces, 
fi is the Poisson's ratio, E is the Young's modulus and parameter a 
is related to the Poisson's ratio by formula a = j^r- The body Q is 



assumed to be elastically isotropic. The equation (1.1.4) is solved for 
the known distribution of body forces in Q and the boundary conditions 
at T defined for displacements (Dirichlet conditions) 

U\ r = $i (1-1.5) 

or for traction forces (Neumann conditions) , 

(a, v) = $ 2 . (1.1.6) 

Here a is a (0, 2)-tensor (called the stress tensor), whose components 
are related to the components of the displacement gradient through 
Hooke's law: 

Oij(U) — (a — 1) 5ij div U + dUi/dxj + dUj/dxi, 

i,j = l,...,n, where dij is the Kronecker symbol. Note that the 
Hooke's law is given above in dimensionless form corresponding to the 
unit value of the shear modulus. The scalar product (a, v) is a vector 
in M. n with the components 



E 



i(U) Uj, i = l,...,n. 



At any part of V the boundary condition can be specified for the dis- 
placement, or for the traction force, or for a linear combination between 
displacements and traction forces. As in the direct electrostatic prob- 
lem, only one boundary condition can be assigned at any part of V. 
The attempt to define simultaneously two different types of boundary 
conditions at the same part of V (i.e. to impose the Cauchy conditions 
|MoFet chapter 6] corresponding to the overconstraining of the system) 
may lead to the loss of the solution. 

The properties of the direct electrostatic and elasticity problems have 
been studied intensively for almost two centuries. It is well-known that 
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problems (1.1.1) and (1.1.4) have unique solutions under the Dirich- 
let boundary conditions ( |1.1.2 ) and (1.1.5), respectively. For Neu- 
mann boundary conditions, solutions exist under additional assump- 
tions J n fdx = fpfads and J n Fdx = f r $2<is, respectively, and are 
unique up to additive constants (see pa l ITG] ). 

There are various analytical and numerical methods to find solutions 
of the boundary value problems for Poisson and Navier equations. Most 
of the numerical methods developed for these problems are based on 
finite difference approximations [H]] IMGj [5a] ISt] , finite element analysis 
|Ba| ISal ICS] and Fourier transform [Du, KhJ.The finite element method 
has become a dominant approach to solving the elasticity problems, 
with the exception of the microelasticity analysis for strain interactions 
in microstructures, where the Fourier transform is still used intensively. 
All major numerical techniques are still used for the electrostatic (or 
magnetostatic) and electromagnetic problems. 

1.2. Inverse problems for Poisson and Navier equations. In the 

present paper we consider the following inverse electrostatic and elas- 
ticity problems: the charge distributions or the internal body forces 
are not known, and the states of the system, i.e., the functions / or J 7 , 



should be determined from the boundary conditions (1.1.2 )— ( 1.1.3 ) or 



respectively, ( 1.1.5 )-( 1.1.6). Note that it is assumed that both Diri ch- 



ief and Neumann data are obtained from the boundary measurements. 
Such problems have attracted much interest in the recent years among 



physicists and engineers (see section 1.4 and references therein). 

The problems described above can be viewed as examples of inverse 
problems of potential theory |Is3j . These problems are different from 
the Calderon's inverse conductivity and elasticity problems, for which 
the coefficients of the left-hand sides of the equations, rather than the 
right-hand sides, are unknown and have to be determined from the 
boundary data (see, for example, [Cal ETJ1 EQ IAMR1 Hs2] ). 

The inverse electrostatic problem formulated above is closely related 
to the inverse gravimetry problerrj^] that has important applications to 
geophysics and has been intensively studied for many years (see, for 
instance, |Isl| IIs2l IMiFoj and references therein). 

In the present subsection we collect some general results on unique- 
ness of solutions of inverse problems for Poisson and Navier equations. 
Essentially, they are well-known (see, for example, |BSBj ). We present 



their proofs in subsection 34^ for the sake of completeness. 

Let, as before, Q C R n be a Euclidean domain with piecewise smooth 
boundary T. Consider the following overdetermined boundary value 



1 We thank Leonid Polterovich for bringing this link to our attention. 
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problem for the Poisson equation: 

Au = /, x en, (1.2.1) 

u|r = 0, (Vw,z/)| r = 0. 

Let H(Q) be the space of harmonic functions on Q. Denote by Z(Q) 
its orthogonal complement in L 2 (Q). We have the following 



Theorem 1.2.2. A nonzero solution of problem (1.2.1) exists if and 
onlyiffeZ{Q). 

Let V C L 2 {VL) be a linear subspace. We say that the inverse elec- 
trostatics problem possesses a uniqueness property for charge distribu- 
tions in V if for any two solutions u and w of the Poisson equations 
Am = / and Aw = g in Q with f,g E V, the equalities u\r = w\r 
and (Vu,v)\r = (Vw,u)\ r imply / = g. Since V is a linear subspace 
of L 2 (Q) and the Poisson equation is also linear, this is equivalent to 



saying that for any nonzero f E V, problem (1.2.1) does not have a 



solution. Therefore, Theorem |1. 2. 2 implies the following 



Corollary 1.2.3. The inverse electrostatics problem possesses a unique- 
ness property for charge distributions in a linear subspace V(Q) C 
L 2 {Vt) if and only ifV(Q) n = 0. 



Remark 1.2.4. It follows immediately from Corollary 1.2.3 that the in- 
verse electrostatics problem possesses a uniqueness property if V(Q) C 
H(Q). For instance, this is true if V(O) is the space of linear functions 
on Q. 

Similar results hold for the inverse elasticity problem. Consider an 
overdetermined problem for the Navier equation: 

AU + a grad div U — F, x E Q, (1.2.5) 
U\ r = 0, (a,u)\ r = 0. 

Let 

£ = A + a grad div 

be the Navier operator acting on vector- valued functions U : Q — > M n . 
Denote by WiVl) the kernel of £ (i.e., the analogue of harmonic func- 
tions for the Navier operator) and by Z(Q) its orthogonal complement 
inL 2 (fi,M"). 



Theorem 1.2.6. A nonzero solution of problem (1.2.5) exists if and 
onlyiffEZiSl). 
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Let V(Q) C L 2 (Q,W n ) be a linear subspace. We say that the inverse 
elasticity problem possesses a uniqueness property for internal stress 
distributions in V if for any two solutions U and W of the Navier 
equations CJJ = F and CW = G in fl with F, G G V, the equalities 
U\r = W\r and (<tu,v)\t = (o~w,v)\r imply F = G. Here ou and 
cr^ denote the stress tensors associated with U and W, respectively. 



Since the Navier equation and the space V(f2) are linear, Theorem 1.2.6 
immediately implies 

Corollary 1.2.7. The inverse elasticity problem possesses a unique- 
ness property for internal stress distributions in a linear subspace V(f2) C 
L 2 (Q,R n ) if and only ifV(ft)nZ(ft) = 0. 

1.3. Checkered distributions. In practical applications one can of- 
ten assume that the distributions of electric charge as well as of the 
internal stress have a certain structure. For example, the geometry of 
the charge density distribution in the electronic component can be dic- 
tated by the structure of the component. Such structures often consist 
of elements with rectangular shape (see, for instance, |LK[ IXCS] ). Let 
us also note that similar structures appear naturally in geophysics |Tsj . 
This motivates the following definition. 

We say that the set IT C M. n is a boxif II = [a 1; 6 X ) x • • • x [a n , b n ), < 
bi, i = 1, . . . , n. Denote by V C (H) C L 2 (I1) a linear subspace generated 
by the characteristic functions of all boxes contained in II. Elements 
of K(n) are called checkered functions. Equivalently, a function / 6 
L 2 (II) is checkered if II can be represented as a finite union of disjoint 
boxes, II = 111 U ■ ■ ■ U IIjv, such that f\ui = const, i — 1, . . . N (such a 
representation is clearly not unique). Note that the subspace V C (JV) is 
dense in L 2 (II). 

Theorem 1.3.1. Let u and w be solutions of the Poisson equations 
Au = f and Aw = g in the interior of the box II C M" with f,g& 
V C (U). Ifu\ 9 u = w\au and (Vw, u)\ v = (Vw, v)\dn, then f = g. 

In other words, the inverse electrostatics problem on II possesses a 
uniqueness property for electric charge distributions given by checkered 
functions. 

Remark 1.3.2. In the context of the inverse gravimetry problem, the 



right-hand side of equation (1.1.1) should be understood as the mass 



density and the function u as the gravitational potential. Therefore, 



Theorem 1.3.1 can be reformulated as follows: the inverse gravimetry 
problem possesses a uniqueness property for mass distributions given 
by checkered functions. To our knowledge, distributions of this type 
have not been previously studied for the inverse gravimetry problem. 
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Uniqueness results for other types of mass distributions could be found 
in |Isl4 Corollary 4.2.3] and |Is3| Theorem 2.1]. 

An analogue of Theorem |1.3.1| holds also for the inverse elasticity 
problem. Denote by V C (I1) C L 2 (U, W 1 ) a linear subspace generated by 
functions F = (/i, f 2 , . . . , f n ), where f t E V C (IL), i = 1, . . . , n. 

Theorem 1.3.3. Let U and W be solutions of the Navier equations 
AU + agraddivf/ = F and AW + agraddivW = G in the interior 
of the box n C E n with F,G e V C (II). Suppose that U\ 9n = W\ 9n 
and (ay, v) |an = (o~w,v)\qu, where ay and aw are the stress tensors 
associated with U and W , respectively. Then F = G. 

In other words, the inverse elasticity problem on II possesses a unique- 
ness property for internal stress distributions with components given 
by checkered functions. 



Theorems 1.3.1 and 1.3.3 are proved using Corollaries 1.2.3 and 1.2.7 
see section [2j 

1.4. Discussion. The interest in the inverse problems considered in 
the present paper arises from a number of practical applications. For 
example, in microelectronics, the observation of the internal voltage 
distribution in a device can be very important for the testing and di- 
agnostic of devices under development |LK| [BSD], 

The inverse elasticity problem naturally appears in the analysis of 
residual stresses |Wil} |Wi2j. These stresses are produced in the ma- 
terials as a result of non-uniform deformation during forming, heat 
treatment and welding processes. The effect of a residual stress field 
is similar to the effect of an internal force distribution, and one can 
be converted into the other. Modern experimental methods, such as 
Scanning Probe Microscopy |KBS[ IGAT] , can be used to obtain data 
on the electric potential and electric field at the surfaces of the compo- 
nent [PrJ . Digital image correlation |CRS] can be applied to study the 
displacement distribution. These methods allow to obtain the Cauchy 
boundary conditions for electrostatic or elastic problems corresponding 
to real objects or components with high accuracy and fine resolution. 
The important question is to which extent such information can be used 
to find the charge (and the potential) or the internal force distributions 
inside the body, and whether the corresponding inverse problems have 
unique solutions. For simple distributions of internal charges or resid- 
ual stresses the inverse problems can be solved easily (for example, 
for a 2-D distribution of charges in a thin layer or 1-D distribution of 
residual stresses with a single significant stress component). However, 
for general charge and internal stress distributions the issue becomes 
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quite difficult. While Theorems 1.3.1 and 1.3.3 give a complete math- 
ematical solution of the inverse electrostatics and elasticity problems 
for checkered distributions, from the viewpoint of practical applications 



these results are far from satisfactory, see section 2.4 



1.5. Non-uniqueness of solutions: an example. One may ask 



whether the analogues of Theorems 1.3.1 and 1.3.3 hold for other, non- 
checkered, electric charge and internal stress distributions. Below we 
provide an example of a natural class of distributions for which the 
solutions of the inverse problems are not unique. Similar examples are 
well-known for the inverse gravimetry problem (see |Is3j ). 

Let S = S(ri,r 2 ) be a spherical layer centered at the origin, that 
is S(ri,rz) = {t\ < \x\ < r 2 } for some r 2 > r\ > 0. Denote 
by V„(S) C L 2 (S) the linear subspace generated by characteristic 
functions of spherical layers centered at the origin. In other words, 
/ G V a (S) if and only if there exists a decomposition of S into a dis- 
joint union of spherical layers S — S± U • • • U Sn, such that f\si = const, 
i — 1, . . . N. We also denote by V a (S) C L 2 (S, R n ) the linear subspace 
of vector functions whose components belong to V a (S). 

Theorem 1.5.1. Let 5cl™ be a spherical layer. Then 
(i) V a {S) n Z(S) ^ {0} and (ii) V*(S) H Z{S) ^ {0}. 



Theorem 1.5.1 is proved in subsection |3.2| Together with Corollaries 



1.2.3 and 1.2.7 it immediately implies 



Corollary 1.5.2. The solutions of the inverse electrostatics and elas- 
ticity problems are not unique in V a (S) and V (T (5') ; respectively. 

1.6. Plan of the paper. Section [2] is devoted to the proof of Theo- 

In subsection |2.1| an auxiliary discretization of 



rems 



1.3.1 and 1.3.3 



the checkered functions is constructed. In subsection 12.21 we introduce 
a family of harmonic functions given by complex exponentials, that are 
used to show that there are no nonzero checkered functions orthogonal 
to the space of harmonic functions. Theorem 1.3.1 then follows from 



Theorem 1.2.2 In subsection 2.3 the above arguments are modified 



in order to prove Theorem 1.3.3 Theorems 1.2.2 and 1.2.6 as well as 



Theorem 1.5.1 are proved in section |3j 



2. Inverse problems for checkered distributions 



We 

which is the most interest- 
ing case for applications. A similar argument works in any dimension 
n > 2. 



The goal of this section is to prove Theorems 1.3.1 and 1.3.3 
present the proofs in three dimensions 
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2.1. Discretization of checkered functions. Let f2 be a box as de- 



fined in section 1.3 For any / G V c (£l), let us construct a function 
/ supported on a finite number of points. Consider an arbitrary box 
IT = [ft" ft+) x [ft", ft+) x [ft-, ft+) c M 3 . Set 

T (Xn) = <Tl(T2(73l (Pi 1 ,^ 2 ,l3s 3 ) (2.1.1) 

Here l^,^) is a function that takes value 1 at the point (x,y,z) and 
vanishes elsewhere. The function T(xu) is supported on the vertices 
of II and takes values ±1 at each vertex. The map T can be then 
extended by linearity to the whole space V^(O) . 

Given a function / G V C (Q), set / = T(/). Denote by V C (Q) the 
space of functions supported on finite subsets of Q. 

Proposition 2.1.2. The map T : V C (Q) — > V C (Q) is injective. More- 
over, there exists a constructive procedure to recover f G V C (Q) from 
the function T(/) = /. 



To prove Proposition 2.1.2 we need an auxiliary lemma below. 

Let {(xi, yi, z{j}f =1 be the collection of vertices of all the boxes ap- 
pearing in some representation of / as a linear combination of charac- 
teristic functions of boxes. We say that a point (a, b, c) G £1 is a node 
of the function / if a = Xi,b — yj,c = for some 1 < i,j, k < N. 
A node v is interesting if f(v) ^ 0. We also call a node v artificial if 
there exists a neighborhood of v in which / does not change its value 
across a plane passing through v and parallel to one of the coordinate 
planes. It is easy to check that all artificial nodes are not interesting 
(and, therefore, artificial nodes can not be determined from /), but the 
converse is not necessarily true. 

Example 2.1.3. Let O be a cube with side 2 centered at the origin 
v = (0,0,0). Let / be a restriction to Q of a function which is identi- 
cally equal to 1 in the positive and the negative octants, and vanishes 
elsewhere. Then v is not an artificial node, but at the same time 



f(v) = by (2.1.1) and, hence, v is not interesting. 



Remark 2.1 A. One could view the difference between artificial and non- 
artificial nodes as follows. Let us colour Q in such a way that points 
x,y G f2 have the same colour if and only if f(x) = f(y). Then Q 
can be represented as a disjoint union of sets Q = Uj =l Qj, such that 
all points in Qj,j = 1, . . . J, have the same colour, and the points in 
Qi and Qk, i 7^ k have different colours. Each Qj is a not necessarily 
connected union of boxes. A node is not artificial if it is a vertex of 
one of the sets Qj, and artificial otherwise. 
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Let supp / = {(pi, qi, be the set of interesting nodes. We say 

that a point (a, b, c) G Q is a marked node if a = p t ,b = qj,c = Sk for 
some 1 < i,j, k < M. 

Note that the properties of being a marked node or an interesting 
node do not depend on the choice of the representation of /. 

Lemma 2.1.5. The set of all marked nodes contains the set of all 

non-artificial nodes. 

Proof. Without loss of generality, suppose that the node (0,0,0) is not 
marked. This means that among interesting nodes there are either no 
points with x = 0, or with y — 0, or with z = 0. In each case, the 
corresponding plane (say, x = 0) does not contain interesting nodes. 
Let us show that the function / does not change its value across this 
plane. This would mean that all nodes contained in this plane are 
artificial, including (0,0,0). 

Consider a decomposition of Q into boxes, such that the set of all 
their vertices coincides with the set of all nodes of / (this could be 
achieved by constructing planes through each node parallel to the co- 
ordinate planes). It follows from the definition of a node that / is 
constant on each of these boxes. Take one of the corner nodes belong- 
ing to the plane x = (i.e. a node lying on one of the edges of fl). 
At each such node at most two boxes meet. Therefore, if this node is 
not interesting, the values of / at the boxes adjacent to it are equal 



and hence the node is artificial. Note that by formula (2.1.1 ), the total 
contribution of these two boxes to the value of / at any other node 
lying on the plane x = is zero. Let us throw away these two boxes 
and pick another node where at most two of the remaining boxes meet. 
Again, the value of / at this node is zero and hence the values of / at 
the boxes adjacent to it coincide. Therefore, this node is also artificial. 
We repeat the procedure until all boxes adjacent to the plane x = are 
thrown away. At each step we get artificial nodes only. This completes 
the proof of the lemma. □ 



Let us now prove Proposition 2.1.2 The proof is based on a similar 



inductive argument as above. We start at a corner box, on which 



formula (2.1.1) allows us to reconstruct in an unambiguous way the 
value of / from the value of / on the corresponding corner vertex. We 
remove that box, move to an adjacent one and repeat the procedure. 



A similar approach will be used again in the proof of Proposition 2.2.5 



Proof. As follows from Lemma 2.1.5 knowing supp / allows us to con 



struct a decomposition of Q, into boxes, whose vertices include all non- 
artificial nodes. We know the values of / at each vertex of these boxes. 
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Let us now reconstruct the value of / at each of the boxes using the 
following inductive procedure. Start with a vertex that is also a vertex 
of Q, and take the box that contains it (there is a unique box with this 
property). Since there are no other boxes containing this vertex, by 



(2.1.1), the value of / at this vertex determines the value of / at the 
box. We subtract the contribution of this box to /, throw away this 
box and take one of the new corner vertices, at which at most two of 
the remaining boxes meet. At each step of this procedure we determine 
the value of / on the corner box, and reduce the number of boxes by 
one. Since the number of boxes is finite, eventually we will determine 
the value of / on each box. □ 

2.2. Exponential functions. Let 

e = e (x) = e(a, 0, x) = e «( >*)+M*,*) 

be a function of the variable x G M 3 , depending on the parameters 
^ a G R, G M 3 , * G M 3 , such that (6, *) = 0, |8| = |*| = 1. It is 
easy to check that e(x) G H(M. 3 ). 

We say that a pair of vectors (G, \J r ) is admissible if the plane it 
generates is not orthogonal to any of the coordinate axes. Set 

P f {a, 0, *) = (/, e(a, 0, x)) = [ f( x ) e <^)+^^) dx. (2.2.1) 

Lemma 2.2.2. Let {vj} be the set of interesting nodes of f G V C (Q). 
Then, for any a^O and any admissible pair (0, \I>) we have: 

P f (a, Q,^) = CJ2 /K) <a, 0, v 3 ), (2.2.3) 

where 

C = C(a, 0, vl/) = 1 f[ (2.2.4) 

OL 0/ -\- t^i 

1=1 1 1 

Note that the constant C is well-defined for any admissible pair 
(©,*)■ 



Proof. The result follows from (2.1.1) by a direct computation of the 



triple integral (2.2.1). □ 



Since any function e(x) is harmonic, the right-hand side in (2.2.3) 



can be computed using the boundary data 0i, 2 of problem (1.2.1 ) by 
Green's formula: 

P f (a, 0, *) = (e{x)4> 2 - ^ <P^j ds 
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Proposition 2.2.5. Knowing the value of Pf(a,Q,^f) for any a ^ 
and any admissible pair (0,\l/) ; one can reconstruct the function f. 

Proof. Let K be the convex hull of supp /. It is easy to see that K is 
a convex polyhedron; let Wj be its vertices. Then, for any 6gR 3 and 
any \1/ chosen in such a way that the pair (0, \I>) is admissible, we have: 

log \P f (a, 0,^)1 

limsup — - — = max(u>j, 0) = max(y, 0), (2.2.6) 

a j yeK 



where the first equality follows from (2.2.3) and the second one from a 
well-known fact that the maximum of a linear functional on a convex 
polyhedron is attained at one of the vertices. At the same time, a 
convex set K can be represented as the intersection of its supporting 
half-spaces: 

K = n 0e iR3{x g M 3 I (x, 0) < max(y, 0)}. (2.2.7) 

y£K 



Therefore, by (2.2.6), we can recover K and, in particular, all its ver- 
tices Wj, j = 1, . . . , N. 

In order to recover the values f(vjj) we use the following procedure. 
Let 0(j) be an external unit normal vector to a plane passing through 
Wj and not intersecting the convex set K (external means here that 
0(j) points to the half-space not containing K). One can easily check 
that in this case 

(Q(j),Wj) > (0(j),O (2.2.8) 
for all k^j,k = l,...,N. 

Choose ^(j) in such a way that the pair (O(j), is admissible. 

We have: 

This allows us to determine the values of / at all vertices Wj, j = 



1,...,N. Using Lemma 2.2.2 we can subtract the contributions of 
these nodes from Pf(a, 0, \1/) and repeat the procedure. Since at each 
step the number of nodes decreases, the number of steps will be finite 
and at the end we will recover all elements of supp / and the values of 
/ at each of these points. □ 



Combining Propositions 2.2.5 and 2.1.2 we obtain the proof of The- 
orem 11.3.11 

Remark 2.2.10. The results of this section generalize in a straightfor- 
ward way to any dimension n > 2. Note that in dimension n = 2 the 
admissibility assumption can be omitted, because for any orthogonal 
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nonzero vectors G, \l/ G K , the denominators in (2.2.4) are automati- 
cally nonzero. 



2.3. Proof of Theorem |1.3.3[ Let us indicate how the proof of The- 
orem 1 1 . 3. 1 can be modified in order to prove Theorem 1.3.3 Let 
F = (A, /a, /a) e Vc(fl) and let F = {fx, f 2 , h) be its discretization 
in the sense of section 2.1 We say that Vj G f2 is a node of F if it 
is a node of one of the functions fi, i = 1,2,3. As before, consider a 
harmonic function 

e := e{a,G,^;x) = e «(©-*)+«0M . 

It is easy to check that 

curl(e(a, 9, x), 0, 0) = (0, a(9 3 + i* 3 )e, -a(0 2 + z* 2 )e) G 

This follows from the fact that e is harmonic and that divcurl = 0. 
Similarly, curl(0, e(a, 9, x), 0) G ft(fi). 
Set 

7^(ct, 6, *) = (F, curl(e(a, 9, x), 0, 0)), 

7>f (a, 9, *) = (F, curl(0, e(a, 0, x), 0)) 

(now (•,•) means the natural inner product in L 2 (fi,lR 3 )). Theorem 
1.3.3 can be now deduced from Proposition |2.1.2 and the following 



analogue of Proposition |2.2.5| 

Proposition 2.3.1. Knowing the value of V (a, O, I = 1,2, for 
any and any admissible (in the sense of section \2jfy pair (9, ty) 

one can reconstruct the function F. 



Proof. Similarly to Lemma 2.2.2, we have: 
Vl{a, 9, viz) = 1 fl £ e(v,)(f 2 (v 3 )(Q 3 + i* 3 )- 



l=i 



/3(«i)(e 2 + i* 2 ))- (2-3.2) 

When choosing the unit vector we will make sure that, apart from 
the admissibility condition, the following condition is satisfied: 

3 ip 2 - 2 ip 3 ^ 0. 
This condition guarantees that if 

/ 2 K)(© 3 + i* 3 ) - / 3 (^)(© 2 + m) = o, 

this automatically implies / 2 (fj) = fs{vj) = 0, and so no term in the 
sum (2.3.2) may "accidentally" vanish. Therefore, the contribution of 
each interesting node will be taken into account. Arguing in the same 
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way as in the proof of Proposition 2.2.5 we can recover the convex hull 
of supp f 2 U supp / 3 . 

Taking V F (a,Q,^f) instead of V F (a,Q,ty) in the argument above 
we recover the convex hull of supp fi Usupp / 3 . Taking a union of these 
two sets, we recover the convex hull of supp F. 

Let Wj be a vertex of the convex hull of supp F. Choose a unit vector 



0(j) satisfying (2.2.8) as in the proof of Proposition 2.2.5 Consider 
two admissible pairs (0(j), ^ fl (j)) and (0j, \l/ 2 (j)). Using (2.3.2) we 
can calculate 

/ 2 K-)(e 3 (jl + ^(i)) - /3W(e 2 (i) + ^(j')), k = 1,2. 

We obtain a system of two linear equations on fziwj) and / 2 (Wj). 
Clearly, we can choose the admissible pairs (0(j), and (Qj, ^ 2 {j)) 

in such a way that the determinant of this system is nonzero. Thus, 
we can compute fziwj) and fzfaj). 

Applying the same argument to V F (a, O, we compute f\{wj). 
Therefore, we have computed F(vjj), and this can be done for any 
vertex of the convex hull of suppF. As in the proof of Proposition 



2.2.5, we subtract the contributions of these nodes from V F (a, Q,^f) 
and V F (a, 0, and repeat the argument. The process will stop after 
a finite number of steps because the number of nodes of F is finite, 
and it decreases at each step. This completes the proof of Proposition 
|2~3~T1 and of Theorem UlUl □ 

Remark 2.3.3. Instead of using the curl in the proof of Proposition 



2.3.1, we could take grade(a, 6, x). Clearly, 

grade(a,^,^;s) G H(£l). 

In this case, for each 0(j) we need to consider three admissible pairs 
(0(i), ^(j; k)), k = 1,2,3, in order to get a system of three linear 
equations on ft(wj), k = 1,2,3. The rest of the proof goes along the 
same lines as above. The advantage of this approach is that it works 
in any dimension, while the curl is defined only in dimension three. 

2.4. Computational challenges. It is convenient to prove Theorems 



1.3.1 and 1.3.3 using the exponential functions introduced in subsection 
2.2 In principle, our proof could be presented as an algorithm that 
allows to reconstruct in a unique way the solutions of the Poisson and 
Navier equations from the corresponding boundary values. However, 
numerical implementation of our approach faces serious computational 
difficulties that we describe below. For simplicity, a 2-dimensional 
example is presented. 
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Figure 1. Geometry of the test problem. Gray areas 
represent charge areas. Thick black line represents the 
convex envelope of charge Thin gray line repre- 

sents the boundary T of the domain Q. 

We have tested the developed algorithm for the solution of the inverse 
electrostatic problem on a rectangle Q C M 2 with boundary T, contain- 
ing two rectangular charged areas (Fig. 1). The boundary conditions 
corresponding to this problem were obtained using the Green function 
method implemented numerically. In other words, we computed a so- 
lution u of the equation Au — f on the whole plane using Green's 
function , and calculated numerically its values as well as the values of 
its normal derivative on T. Here / is the characteristic function of the 
total charged area. 

It was found that if P/(a, 0, \E r ) is calculated directly using the in- 



tegration over f2, then the procedure based on (2.2.6) and (2.2.7) pro- 
duces the convex hull of the charged areas with high accuracy. However, 
when Pf(a, 0, \E r ) is calculated using the integration over the boundary, 



the algorithm based on (2.2.6) produces the convex hull occupying the 
whole Q. Such a drastic difference is the result of small numerical er- 
rors in the boundary conditions determined by the numerical solution 
of the direct problem, and also in the numerical integration over the 
boundary. 

Strong sensitivity to numerical errors arises from a specific nature 
of exponential functions. When a is large, these functions are rapidly 
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increasing in one direction and rapidly oscillating in the orthogonal 
direction. As a result, small errors in boundary conditions are mul- 
tiplied by large factors and the convex hull estimate based on (2.2.7) 
and ( 2.2.9[ ) becomes distorted. For small values of a, the value of 
log \Pf(cx, 0, *f?)\/at obtained using the integration over Q is close to 
the one obtained using the integration over T. However, as a increases 
these two values diverge. 

In order to verify the validity of the computational model, calcula- 
tions were performed for the constant charge density: / = 1 on Q. 
In this case, the precise value of Pf(a, 0, \&) can be computed ana- 
lytically. For all values of a, the analytical equations produced the 
same log \Pf(a, 0, ^)\/a values, no matter if the integration was per- 
formed over Q or its boundary T. However, in the numerical analysis of 
this problem, the area and boundary integration were producing close 
values for small values of a and were diverging for large a. 

In practical problems based on experimental data, boundary condi- 
tions are always obtained with some errors. Therefore, as the above 
analysis show, in order to produce an algorithm for the solution of the 
inverse electrostatics problem that is numerically implementable, one 
needs to modify our approach. One possibility would be to find a set of 
harmonic functions exhibiting good behavior from the numerical view- 
point, which could replace the exponentials in the proof of Theorem 

on 



3. Necessary and sufficient conditions for uniqueness of 

solutions 

3.1. Proofs of uniqueness criteria. The goal of this subsection is 



to prove Theorems 1.2.2 and 1.2.6 Let us start with Theorem 1.2.2 



To prove necessity, suppose that u is a solution of problem (1.2.1) and 
let h e H{Q). Then 



f-h 



Au ■ h 



u- Ah = 0. 



Note that the boundary terms in the integration by parts disappear 
since 4>i — 4>2 = 0. 

To prove sufficiency, denote by Gi(x,y) the Green's function of the 
Dirichlet boundary value problem in Q: 



A x G 1 {x,y) = 8(x - y), G 1 (x,y)\ xegn = 0, 
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and by G 2 {x,y) the Green's function for the corresponding Neumann 
boundary value problem: 

A x G 2 (x,y) = 5{x-y) - (gmd x G 2 (x,y),v)\ XIEdU = 0. 

Note that the integral over Q of the right-hand side of the equation 
above is zero, which is necessary for the existence of a solution of the 
Neumann problem with zero boundary conditions. 
It follows from the definitions of G\ and G 2 that 



{G 1 -G 2 )(x,y) + 



2im 



is a harmonic function of x. Therefore, the assumption / G Z implies 



u (y) ■= / f(x)G 1 (x,y)dx = 
for all y G f2. Note that the term 

/(*) 



f(x)G 2 {x,y) dx - / f(x 



XT 



2\n\ 



dx 



x\ 
2\tt\ 



dx 



is constant and hence of no importance for the Neumann boundary 
value problem. Let us also remark that j Q f(x)dx = since / G Z(Q). 
It is easy to check that, by the properties of Gi and G 2 , the function u 



constructed above is a solution of problem (1.2.1). This completes the 
proof of Proposition |1.2.2| □ 



The proof of Theorem 1.2.6 is completely analogous. We note that 
the existence of Green's functions for the Navier equation with either 
Dirichlet or Neumann boundary conditions is well-known — see, for 
instance [Sol section 7.12]. 



3.2. Proof of Theorem 1.5.1 Proof of (i). Let A, B, S be three 
spherical layers centered at the origin such that A U B = S and 
A R B = 0. Consider the following linear combination of charac- 
teristic functions of the sets A and B: 

u = Vo\{B) xa-Vo\(A) X B. 

By the mean value theorem for harmonic functions we immediately 
have u G Z(S), and this completes the proof of part (i) of the propo- 
sition. 



Proof of (ii) In order to prove the second part of the proposition we 
note that a function U lying in the kernel of the Navier operator C 
is biharmonic. Indeed, let C(U) = AU + agrad divU = 0. Taking 
the divergence on both sides we get div(AU) = 0. Here we took into 
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account that div grad = A and that the Laplacian commutes with the 
divergence. At the same time, applying the Laplacian to C{U) we get 

A 2 U + a A grad div U = 0. 

But since div(A[7) = 0, the second term vanishes, because 

A grad div U = grad div AU. 

Hence, A 2 U = and U is biharmonic. 

It is well-known that a real valued biharmonic function f(x) satisfies 
the following mean- value property (see, for example, [EKJ): 

n+2 

f(y)dy = u n r n f(x) + ^—Af(x), (3.2.1) 

B(x,r) l \ n + l ) 

where B(x,r) is a ball of radius r centered at x G 1R 71 and u n is the 
volume of the unit ball in IR n . 

Consider now a ball B = B(r 3 ) centered at the origin. Let us repre- 
sent it as a union of three sets B{ri) U S(r 1 ,r 2 ) U 5 , (r 2 ,r 3 ), r 3 > r 2 > 
?"i > 0. Let F a £ be a piecewise constant vector- valued function taking 
the values a and b on S(r 2 , r 3 ) and S(ri, r 2 ), respectively, and the value 
one on B(ri). Clearly, F a ^ G for all a, 6 e E. Let us show that for 
any triple < r% < r 2 < r 3 there exists a choice of parameters a, b such 
that F 0) 6 G ZiB). It can be deduced from the mean- value formula 



(3.2.1), that the inclusion F a ^ G 2(B) holds if the following system of 



equations is satisfied: 

a(^ 3 -r 2 )+^2-<)+< = 0, a(r 3 1+2 -r" +2 )+6(r" +2 -< +2 )+ri l+2 = 0. 

One may check that the determinant of this system does not vanish 
for r 3 > r 2 > r\ > 0. Indeed, if r 3 ^ r 2 , the only positive roots of the 
determinant considered as a polynomial in r\ are r\ = r 2 and r\ = r 3 
(there are no other positive roots because, as can be easily verified, 
the derivative with respect to r\ has only one positive root). Hence, 
there always exists a unique solution a, b of the system above, and the 
corresponding function F^ b G 2(B). 



This completes the proof of Theorem 1.5.1 



Remark 3.2.2. The proof of part (i) of Theorem 1.5.1 shows that the 
intersection V a (S) D 2(S) is in fact quite large. Indeed, it is easy to 
show that for any partition of S into k concentric spherical layers, there 
exists a (k — l)-dimensional linear subspace of functions in V a R 2(S). 

The proof of part (ii) goes through without changes if instead of 
a ball B(r 3 ) one takes a spherical layer S = S(r 3 ,r ) = S(r ,ri) U 
S(r±, r 2 ) U S(r 2 , r 3 ). It follows from the proof that for any partition of 
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S into k concentric spherical layers, there exists a (k — 2)-dimensional 
linear space of functions in V a fl Z{S). 
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